function [ Fgrid ] = Bivariate_4pointRefinement( grid , omega )
%-----------------------------------------------------------------
% Input:
% grid(m,n) - a grid over Z^2
% Output:
% Fgrid(2*m-5,2*n-5)
%-----------------------------------------------------------------
% ABSTRACT
% Tensor product. the middle point is insert in the order y->x of
% calculation
%-----------------------------------------------------------------
% NIr Sharon, 07-05-12
%-----------------------------------------------------------------

m = size(grid,1);
n = size(grid,2);

Fgrid = zeros(2*m-5,2*n-5);

Firstline = zeros(1,2*n-5);
Lastline  = zeros(1,2*n-5);

for x=2:(m-2)
    for y=2:(n-2)
        Fgrid(2*x-3,2*y-3) = grid(x,y);
        Fgrid(2*x-2,2*y-3) = FourPointAvarage([grid(x-1,y) , grid(x,y) , grid(x+1,y) , grid(x+2,y)],omega);
        Fgrid(2*x-3,2*y-2) = FourPointAvarage([grid(x,y-1) , grid(x,y) , grid(x,y+1) , grid(x,y+2)],omega);
    end
end

% Bounderies handeling ----
for y=2:(n-2)
    Firstline(2*y-3) = grid(1,y);
    Firstline(2*y-2) = FourPointAvarage([grid(1,y-1) , grid(1,y) , grid(1,y+1) , grid(1,y+2)],omega);
    Lastline(2*y-3) = grid(m,y);
    Lastline(2*y-2) = FourPointAvarage([grid(m,y-1) , grid(m,y) , grid(m,y+1) , grid(m,y+2)],omega);
end    
%--------------------------

for x=2:(m-2)
        Fgrid(2*x-3,2*n-5) = grid(x,n-1);
        Fgrid(2*x-2,2*n-5) = FourPointAvarage([grid(x-1,n-1) , grid(x,n-1) , grid(x+1,n-1) , grid(x+2,n-1)],omega);
end
Lastline(2*n-5) = grid(m,n);

for y=2:(n-2)
        Fgrid(2*m-5,2*y-3) = grid(m-1,y);
        Fgrid(2*m-5,2*y-2) = FourPointAvarage([grid(m-1,y-1) , grid(m-1,y) , grid(m-1,y+1) , grid(m-1,y+2)],omega);

end

Fgrid(2*m-5,2*n-5) = grid(m-1,n-1);

%Fy = Fgrid;


for x=3:(m-3)
    for y=2:(n-2)
        Fgrid(2*x-2,2*y-2) = FourPointAvarage([Fgrid(2*x-5,2*y-2) , Fgrid(2*x-3,2*y-2) , Fgrid(2*x-1,2*y-2) , Fgrid(2*x+1,2*y-2)],omega);
    end
end

% Bounderies handeling-----
x=2;
for y=2:(n-2)
    Fgrid(2*x-2,2*y-2) = FourPointAvarage([Firstline(2*y-2) , Fgrid(2*x-3,2*y-2) , Fgrid(2*x-1,2*y-2) , Fgrid(2*x+1,2*y-2)],omega);
end


x=m-2;
for y=2:(n-2)
    Fgrid(2*x-2,2*y-2) = FourPointAvarage([Fgrid(2*x-5,2*y-2) , Fgrid(2*x-3,2*y-2) , Fgrid(2*x-1,2*y-2) , Lastline(2*y-2)],omega);
end
%--------------------------


% for x=3:(m-2)
%     for y=3:(n-3)
%         Fy(2*x-2,2*y-2) = FourPointAvarage([Fy(2*x-2,2*y-5) , Fy(2*x-2,2*y-3) , Fy(2*x-2,2*y-1) , Fy(2*x-2,2*y+1)],omega);
%     end
% end
       

end

